###################################################
## Demographic Analysis
## SCRIPT 1: Run simulations for population 
## tolerance tests.
## 
## Script to run SMC at various population 
## tolerances under Census and both DAS versions.
##
## 5/17/21
###################################################

library(ppmf)
library(geomander)
library(tidyverse)
library(sf)
library(redist)

la_all <- readRDS("../../data/LA/la.Rds") %>%
  rename(sld_upper = sld_up,
         sld_lower = sld_low) %>%
  filter(sld_lower != 106)

la_map <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = pop,
    pop_tol = 0.05,
    existing_plan = sld_upper
  )

la_map4 <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = v4_pop,
    pop_tol = 0.05,
    existing_plan = sld_upper
  )

la_map12 <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = v12_pop,
    pop_tol = 0.05,
    existing_plan = sld_upper
  )

la_map19 <- la_all %>%
  st_set_crs(26915) %>%
  st_transform(crs = 4269) %>%
  redist_map(
    total_pop = v19_pop,
    pop_tol = 0.05,
    existing_plan = sld_upper
  )

# each population tolerance level
vals <- c(0.001, 0.005, 0.01, 0.05, 0.1, 0.2, 0.3, 0.5)

for (i in 1:length(vals)) {
  
  cat(paste("Run", i, "\n"))
  
  la_map <- set_pop_tol(la_map, vals[i])
  la_map4 <- set_pop_tol(la_map4, vals[i])
  la_map12 <- set_pop_tol(la_map12, vals[i])
  la_map19 <- set_pop_tol(la_map19, vals[i])
  
  plans_smc <- redist_smc(la_map,nsims = 2500)
  
  filename <- paste("../../data/LA/sim/parity/smc_", vals[i], ".rds", sep = "")
  saveRDS(plans_smc, file = filename)
  rm(plans_smc)
  
  plans4_smc <- redist_smc(la_map4, nsims = 2500)
  
  filename <- paste("../../data/LA/sim/parity/smc4_", vals[i], ".rds", sep = "")
  saveRDS(plans4_smc, file = filename)
  rm(plans4_smc)
  
  plans12_smc <- redist_smc(la_map12,nsims = 2500)
  
  filename <- paste("../../data/LA/sim/parity/smc12_", vals[i], ".rds", sep = "")
  saveRDS(plans12_smc, file = filename)
  rm(plans12_smc)
  
  plans19_smc <- redist_smc(la_map19, nsims = 5000)
  
  filename <- paste("../../data/LA/sim/parity/smc19_", vals[i], ".rds", sep = "")
  saveRDS(plans19_smc, file = filename)
  rm(plans19_smc)
}
